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Three bosonic, spin-polarized atoms in a spherical oscillator potential constitutes the simplest nontrivial Bose- 
Einstein condensate (BEC). The present paper develops the tools needed to understand the nature of the complete 
£\j > J — energy spectrum for this prototype system, assuming a sum of two-body potentials. The resulting spectrum 
| is calculated as a function of the two-body scattering length a sc , which documents the evolution of certain many- 
body levels that evolve from BEC-type to molecular-type as the scattering length is decreased. Implications for the 
behavior of the condensate excited-state spectrum and for condensate formation and decay are elucidated. The energy 
levels evolve smoothly, even through the regime where the number of two-body bound states iVj, increases by 1, and 
p^j ' a sc switches from — oo to +oo. We point out the possibility of suppressing three-body recombination by tuning the 
two-body scattering length to values that are larger than the size of the condensate ground state. Comparisons with 
mean-field treatments are presented. 
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I. INTRODUCTION 

Mean-field methods are commonly used to characterize Bose-Einstein condensates (BECs), i.e., atomic alkali vapors 
in an external trapping potential jlj . These approaches, which account for the interaction between particles through 
a mean-field term in an effective single-particle Hamiltonian, accurately describe much of the dilute condensate 
energetics. However, these shape-independent s-wave scattering length approximations break down for strongly 
interacting gases with large interaction parameters Jl| ^ . 

Very recent experiments by Greiner et al. |^,^] have entered this regime: A quantum phase transition from a 
superfluid to a Mott insulator was observed for ultracold atoms held in a three-dimensional optical lattice. Most 
relevant to our studies here, the atoms are held in tight isotropic lattice sites with an occupancy of 1-3 atoms per 
site. The experiment thus realizes a strongly correlated many-body (and in a sense, few-body) quantum system with 
unprecedented control of parameters. In this paper, we present the quantal J = energy spectrum of three particles 
\ under external isotropic confinement, and compare with results obtained using mean-field theory. 
f^} • The lowest order mean-field approximation [Gross-Pitaevskii (CP) or Hartree-Fock (HF) equation] Jl|J^-[ll|] and 
CnJ \ also higher order approaches such as second order perturbation theory (ljjl^l treat the metastable BEC state of the 
trapped alkali vapor as if it was the "true" ground state of the system. This means that these approaches cannot 
describe recombination processes or atom- loss processes, since coupling to the true (solid-like/liquid-like) ground 
state is neglected completely. Attempts to include additional "empirical" terms in these mean- field equations, which 
account for "reaction processes", exist, and have been able to describe some experimental findings successfully [fl4| [l6| . 
However, the ad hoc nature of some of these empirical terms is less than satisfactory. 

The nonperturbative quantal calculations permit us to circumvent these problems. These calculations come at a 
price, namely, our initial restriction to only N — 3 particles with total angular momentum J = in the trap. For this 
"simple" system, however, we are able to include all channels in the calculation, molecular three-body bound state 
channels, diatom plus atom channels, and metastable gaseous-like channels |T?j] . 

This paper focuses on regimes where standard mean- field treatments break down. For large positive a sc — > +oo, 
the mean-field treatment results in unphysically large diverging energies. In contrast, the three-body energies change 
smoothly in the pole region, where the scattering length a sc changes from — oo to +oo and the number of two-body 
s-wave bound states Ni, increases by 1. Additionally, our three-body study implies that the excitation frequencies are 
well behaved around |a sc | — > oo. 

The mean-field treatment also exhibits limitations in the parameter range around the negative critical scattering 
length a cr , where a cr is defined through the instability criterion (N — l)a cr /a;j = —0.575 derived from the mean-field 
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equation 11 lq]. The lowest gaseous energy level (i.e., the BEC ground state) changes its character from metastable 
to molecular around a cr , while some of the higher lying BEC-like states remain almost completely unaffected by 
criticality. In contrast, the mean-field equation simply has no solution for a sc lower than the critical value. 
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Three-particle studies have been used previously to shed light on the many-body physics of BECs. In 1996, Esry et 
at explored the role of the two-body scattering length in BECs. More recently, three research groups explained 
experiments regarding the three-body recombination rate in many-atom BECs via three-body calculations that include 
channel coupling [[l7|]2C| |23]| . Three-body studies have also been used to test the validity of shape- independent two- 
body interaction potentials p4|] . 

The studies presented here are based on time-independent Schrodinger wave mechanics. The energy levels show 
several avoided crossings when plotted as functions of a sc (see Section HI); nevertheless, dynamical behavior such as 
Landau-Zener-type branching ratios or recombination rates remains beyond the scope of this paper. Studies along 
these lines will be a natural continuation of the static studies presented here. 

Section [n] describes our three-body system, and outlines the formalism for its investigation. Section [II presents 
quantum mechanical three-body energetics, and Section IV compares with results derived from mean-field treatments. 
Finally, Section ^ summarizes our main results. 



II. THE THREE-PARTICLE SYSTEM 

The Schrodinger equation for three mass-m particles in a spherical external potential with trapping frequency 
Vho = Wao/(27r) reads 



^ 2m r * 

i=i 
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*(n,r 2 ,r 3 ) = E*(fi,f 2 ,r 3 ), 



(1) 



where V(rij) denotes a two-body interaction potential and r,j = \fj — fi\. fi is the position vector of atom i relative 
to the center of the trap. Restriction to states of zero total angular momentum (J = 0) reduces this nine-dimensional 
problem to six dimensions. The dimensionality of the problem can be reduced further by a transformation to Jacobi 
coordinates, X cm = (ri + r*2 + ?%)/3, pn = T% — Ti, and (012,3 = r 3 — (f2 + ri)/2 [p4| . This transformation decouples 
the center-of-mass motion, 
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Here, M = 3m is the total mass, and = (k + 3/2)^^0 with k = 0, 1, 2, • • ■ the energy related to the center-of-mass 
motion. The corresponding eigenfunctions 4>kim are simply the three-dimensional isotropic oscillator solutions. Eq. (|^) 
represents the Hamiltonian describing the relative particle motions in terms of three hyperspherical coordinates (R, 
i9, tp) rather than the Jacobi coordinates P12 and p*i2,3 [p4|- p7| . Our definition of the hyperradius R, and the two 
hyperangles d and ip is based on the democratic coordinate system of Whitten-Smith [^8), i.e., fiR 2 — mp 2 2 /2 + 
2mp\ 2 3/3. Owing to the identical boson character of the atoms, the angular range can be restricted to ip e [0,7r/6] 

and $ G [0,71-/4]. p — m/y/3 denotes the reduced mass of the three-body system, while A 2 denotes the "squared 
grand angular momentum" operator ^9|. Here we rescale the wave function ipi^R,!),^) by i? 5 / 2 , i.e. ^{^1,7%, 7%) = 
R~ 5 / 2 ?p(R, (p)(j)kim{X) , hi order to eliminate first derivatives from the hyperradial kinetic energy operator. 

Eq. (||) is solved here using two different, but related approaches: a coupled-adiabatic-channcls approach and a 
strict adiabatic approximation. The coupled-adiabatic-channel calculation can in principle be made exact. The strict 
adiabatic approximation can be viewed as an "incomplete" coupled-adiabatic-channel calculation, truncated to one 
channel only, which neglects all off-diagonal coupling elements. It is important to keep in mind that the strict adiabatic 
approximation cannot describe recombination processes. 

To solve Eq. (||) we expand ij) n (R,'d, ip) into radial wave functions F l/n (R) and a set of complete, orthonormal 
angular channel functions <fr u (R; t?, ip). The basis functions <!>„, which depend parametrically on i?, are solutions of 
the partial differential equation 
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<S> v {R-$^) = U v {R)<S>„{R;d,<p). 



(4) 
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We refer to the eigenvalues U U (R) as adiabatic potential curves. Inserting the expansion tjj n = ^2,F vn § u into the 
Schrodinger equation for the internal motion [Eq. (|)] results in an infinite set of coupled ordinary differential equations, 

/r l -^ 2 ho R 2 + U V {R) + Q VU (R) - E*A F vn (R) = (5) 



2(i dR 2 T 



where P vv i [Q U u'] are angular coupling matrix elements involving the first [second] derivative of with respect to 
R |||. We refer to the solution of Eq. (|) as the coupled-adiabatic-channel solution. The strict adiabatic approxi- 
mation follows immediately by neglecting coupling elements between different channels [= setting the right hand side 
of Eq. (||) to zero] . Another variant of the strict adiabatic approximation additionally neglects the diagonal element 
Qw, and is typically referred to as the hyperspherical Born-Oppenheimer approximation. In the following, we use 
these two variants of the strict adiabatic approximation interchangably, and denote the eigenvalues by E™£ . 

The outlined formalism splits the solution of the Schrodinger equation for the internal motion [Eq. (0)] into two 
steps: i) solution of Eq. (^) (here via a two-dimensional B-spline code), and ii) solution of the one-dimensional coupled 
equations, Eq. (||). The former is numerically more challenging due to its higher dimensionality, and coincides with 
determining the bound and continuum channel functions $„ of three particles interacting via a sum of two-body 
potentials V(r) [the trapping potential does not enter Eq. (Q)]. 

For a vanishing interaction potential, V(r) — 0, the eigenfunctions $„ of Eq. (Q) reduce to Gegenbauer polynomi- 
als [^9| with eigenvalues U U (R) — fi 2 A ^ A "^jy 1 , A = 0, 4, 6, • • •, where the A = 2 state is forbidden due to symmetry 
constraints, and the eigenvalues with A = 12 and A > 16 are doubly-degenerate. Now consider a non-vanishing two- 
body interaction potential V(r). Compared to the V = potential, the non- vanishing V introduces three-body bound 
states at small hyperradii R, accounting for the short range physics of V(r), and also modifies the potential curves at 
large R, reflecting the non-zero s-wave scattering length a sc . Additionally, a non- vanishing V lifts the degeneracy of 
the U V {R) potential curves. 

To illustrate this behavior, Fig. |l| shows the sum of the trapping potential, Vt rap = ^h R 2 /^ with i>ho = 78kHz, 



and the adiabatic potential curves, U V (R), for a two-body model interaction potential (see Sec. Ill) with s-wave 
scattering length a sc — 228a. u. and two two-body s-wave bound states, AT b = 2 [m = m( 87 Rb)]. Length and energy 
are expressed in oscillator units (length unit a] lo = y l /h/(muJho), energy unit TiiOho)- The two lowest potential curves 
with v — and 1 (upper panel, dashed lines) describe three-body bound state physics, and approach the two-body 
s-wave binding energies in the absence of a trapping potential. The higher lying potential curves (solid lines) describe 
metastable BEC physics, and would correspond to continuum states, in the absence of a trapping potential. The 
v = 2 curve (thick solid line) shows a double minimum structure with minima at R ss 0.2a^ o and R « 2.5a^ , which 
are separated by a potential barrier. The lowest metastable state "lives" in the minimum at larger R. The lower 
panel of Fig. [I] shows an enlargement of the gaseous region together with the effective potential curves for a vanishing 
interaction potential (dotted lines). Significant deviations between the solid and the dotted curves (non- vanishing 
and vanishing two-body potential, respectively) are visible. 

As clearly demonstrated in Fig. [l], to properly describe the metastable states, one must solve Eq. (^|) for hyperradii 
R as large as 10a^ o , corresponding roughly to several 10 4 a.u. (depending somewhat on V} lo ). In our calculation, we 
determine the adiabatic potential curves U v and the coupling matrix elements P vv i and Q vv i out to R w 10 — 100a sc 
using a two-dimensional B-spline code, and then extrapolate to lar ger R. For the continuum potential curves we follow 
a second independent approach suggested by Nielsen and Macek [g0| . For large R, they recognise the dependence of 
U u on a sc only, and derive a nearly analytic expression for U u . Note, Nielsen and Macek's formula does not account 
for those potential curves, which simply coincide with the "unperturbed" eigenvalues h 2 ^^J^ 5 ^ 4 > A = 12, 16, 18, • • • 
(see above) . The nearly analytical expressions |2~C| are in excellent agreement with our continuum curves U v calculated 
through numerical solution (and extrapolation) of Eq. 

In passing, we note an additional nice feature of the outlined formalism. Suppose that we wish to study the 
energetics of three particles in an external spherical potential as a function of the trapping frequency Uho- For a given 
interaction potential V(r), it is only necessary to solve the numerically most demanding Eq. (||) once. u>ho then enters 
the simple one-dimensional Eq. (^J) simply through an additive term, which trivially determines the energy spectrum 
versus u>} lo . 
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III. THE QUANTAL ENERGY SPECTRUM 



This section presents the internal energies E mt for three bosonic particles with J = in an external trap for 
a two-body model potential of the form V(r) = e?cosh~ 2 (r/r ). This is convenient because the two-body s-wave 
eigenenergies and eigenfunctions can be determined analytically 31 1. Furthermore, we determined an analytical 
expression for the energy-dependent s-wave scattering length, 

a sc (k) = (6) 



where 



ik 1 



iF\ (ik — s, 1 + s + ik;l + ik 



■ 2i 



2* k 2 F 1 {-s,l + s-l-ik;l) 



(7) 



where k — y mE/h , and s — ~ 5 + 1 + ^mdr^jU . 2F1 denotes the usual hypergeometric function. The zero- 
energy s-wave scattering length a sc is then simply given by a sc = lim^o a sc (k). While ro is fixed at ro = 55a. u. 
throughout our calculations, d is varied to change the number of two-body s-wave bound states Nb and the two-body 
s-wave scattering length a sc [m — m( 87 Rb) throughout the rest of the paper]. We perform three-body calculations 
for 70 different d values, which translates into a range of two-body scattering lengths a sc from — 10 4 a.u. to 10 4 a.u.. 
To keep the coupled-adiabatic-channel calculations tractable we restrict d in this study to values such that Nb = 1 or 
2. 

Consider a vanishing interaction potential, d — (ideal gas) first. For this system, the strict adiabatic approximation 
is exact, and E" lt = 3, 5, 7, • • • fruJho- Introduction of a non-vanishing interaction potential shifts these gaseous energy 
levels according roughly to the value of a sc , lifts the degeneracy, and additionally gives rise to three-body bound 
states. Figure |] shows the internal energies -E™' (solid lines) of the gaseous-like states for Uho = 78kHz calculated 
via the strict adiabatic approximation for various scattering lengths a sc . To accomodate a large range of interaction 
parameters (N —l)a sc /ai lo pH , Fig. |] shows the internal energies as a function of [arctan(a sc /a/ 10 )]/7r. To additionally 
be able to plot results for systems with Nb = 1 and Nb = 2 two-body bound states, respectively, the integer part of 
the abcissa is chosen to be equal to — Nb. This scaling procedure results in abcissa values of [—2.5, —0.5]. We refer to 
the interval [—0.5, —1.5] as the "first cycle" (Nb = 1), and to the interval [—1.5, —2.5] as the "second cycle" (Nb = 2). 
In this representation, a sc = corresponds to [arctan(a sc /a/j )]/7r = —2.0 and —1.0. 

To interpret this complicated energy level structure, concentrate on six different scattering length values labeled by 
A-F in Fig. ||. Figure || shows the sum of the corresponding adiabatic potential curves and the trapping potential, 
U v + Vtrap, where V tra p(R) = 0.5/iw 2 o i? 2 (solid lines): 

Panel A corresponds to a large positive scattering length, a sc — 2138a. u.. The lowest BEC-like curve (v — 2) is 
highlighted as a thick solid curve. The lowest energy level lying within this potential curve, -E^q , has an energy close 
to 5huJho, and the excited levels (E^, n > 0) lie about 2hajho higher each (see Fig. ^). We refer to the energy levels 
"located" within the same potential curve (identical is, but having a different number of hyperradial nodes n) as a 
family. The lowest energy of the next higher BEC-like curve is E%q w Ihio, and again, the excited levels (E^ , n > 0) 
lie about 2hu>ho higher each. To guide the eye, Fig. || displays dotted curves that represent the effective potential 
curves in the absence of any two-body interaction potential. 

Panel B corresponds to a much smaller scattering length than in panel A, a sc — 103a. u.. This panel shows a series of 
crossings at energies around 20hujho, indicating a <i-wave shape resonance. As the scattering length decreases slightly, 
these crossings move to lower energies, until they finally vanish at a sc = 97a. u.. The scattering length a sc = 97a. u. 
corresponds to a well depth d of the two-body potential V at which V first supports a d-wave two-body bound state. 
To calculate energy levels in the presence of a d-wave shape resonance we follow the potential curves that cross 
diabatically. Note that the potential curves approach the V(rij) = curves (i.e., the dotted lines) as a sc approaches 
zero. Accordingly, the energy levels approach those of the ideal inhomogeneous gas. 

Panel C shows the effective potential curves (solid lines) for vanishing scattering length, a sc = (however, non- 
vanishing V), in which case the potential curves are almost identical to the ones for vanishing V (dotted lines). Small 
differences occur in the inner wall region reflecting the existence of molecular states at small R. These molecular 
states at small R are, however, not shown in Fig. g|. 

Panels D, E, and F correspond to negative scattering lengths. The inflection of the thick solid line at R « 0.5a/j O 
in panel D indicates the first signature of a decaying condensate. The interaction parameter for this system is 
(iV— l)a sc /aho = —0.085. In panel E, the barrier of the lowest gaseous-like curve, which separates the BEC minimum 
from the energetically much deeper molecular minimum, has almost vanished; here, (N — l)a sc /a/ lo = —0.356. Panel 
F shows a large negative scattering length case, a sc — — 1970a. u.. The thick solid curve corresponds to a situation 
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in which the gaseous-like state (condensate) is no longer stable. Note the similarity of the potential curves for large 
positive and large negative a sc (panel A and F). This similarity reflects the smooth dependence of £™f as |a sc | goes 
through infinity (see below). 

In addition to the aforementioned issues, four more characteristics are important to the energy level scheme shown 
in Fig. ||: (i) At a sc = 0, corresponding to [arctan(a sc /a/ lo )]/7r = —I and —2, the positive energy levels E l j£ lie very 
close to those of the ideal gas. (ii) The energy level with energy E mt = 15huiho is unchanged over the two cycles 
shown. This energy level originates in the "unperturbed" potential curve h 2 A ^ A ~^j^ 15 ^ 4 > A = 12 (see Section ||), and 
therefore coincides with an ideal gas level. (in) The internal energy levels in the first cycle (Nt, = 1) and in the second 
cycle (Nb = 2) behave very similarly indicating that the energies are, to a good approximation, independent of the 
number of two-body s-wave bound states A^,. (iv) Notice in particular that the energy levels E % J™ change smoothly 
as |a sc | passes through infinity 

Note as well that families of plunging energy levels with molecular character are present, which to simplify our 
presentation are shown only partially in Fig. |[ Consider the lowest gaseous-like energy level of the first cycle, i.e. 
with Nb = 1. As the gaseous-like state becomes unstable while a sc decreases this energy level decreases through zero 
(at [arctan(a sc /a? lo )]/7r « —1-2), which reflects the transition of this BEC state to one of molecular character. The 
other energy levels of this family also change character from gaseous-like to molecular-like as the two-body potential 
well depth d becomes deeper. For clarity, Fig. || shows these energy levels only for [arctan(a sc /a^ )]/7r > —1.75. For 
each additional two-body bound state, a "manifold" of plunging energy levels is present, reflecting the fact that the 
number of three-body states with negative energy increases far more rapidly than those for the two-body system. 

So far we have discussed results for a trap with trapping frequency vi lo — 78kHz. In addition, we perfor med calcu- 
lations for trapping frequencies that are 10 times larger, and also 10 or 100 times smaller. Since a/ lo = yfi/ (too^o), 
the interaction parameter (JV — l)a sc /aho scales as y/uJho- Adopting the same two-body interaction potentials as 
before, a trapping frequency smaller than v^o — 78kHz therefore leads to a narrower [arctan(a sc /a/ lo )]/7T interval 
than shown in Fig. ||. As an example, Fig. ^ shows the internal energies E"^ for three particles in a trap with 
Vho = 780Hz and 2V& = 2 (solid lines). A few of the internal energies at large negative scattering length, corresponding 
to [arctan(a sc /a?i )]/7r rs —2.08, are indicated by diamonds. These data points are not connected by solid lines since 
in this region our grid in a sc is too coarse to allow for interpolation. 

To summarize, our microscopic studies of the energetics of three confined particles as a function of Vho reveal that 
the overall behavior of the energy level scheme (as shown in Fig. ^) remains unchanged as V} lo either decreases or 
increases. As discussed in the previous paragraph, however, Vho determines the range of the interaction parameter. 
Specifically, the behavior of the energy levels around a sc « a cr depends on v^o- Figure [2] indicates a gradual decrease 
of the condensate energy as the scattering length decreases through the region [arctan(a sc /aft )]/7r « —2.1. As v\ w is 
decreased, the separation between the molecular-like region at small R and the gaseous-like region at large R becomes 
more pronounced, and the energy decrease (change from gaseous-like to molecular-like) takes place over a smaller 
range of scattering lengths. In other words, the potential barrier (shown in panel D and E of Fig. [3] for Vho = 78kHz) 
spatially separates the molecular-like and gaseous-like region more for small V] lo than for large Vho- 

The internal energies discussed in this section have been calculated within the strict adiabatic approximation. To 
assess the accuracy of this approach we additionally performed coupled-adiabatic-channel calculations. Figure^ shows 
the energies resulting from the coupled-adiabatic-channel calculation including four channels (pluses) together with 
those resulting from the adiabatic calculation (dotted, solid, and dashed lines) for V} lo = 78kHz. Figure |B| indicates 
excellent agreement between these two sets of calculations for small a sc (left hand side of the figure), however, some 
differences occur at larger a sc (right hand side of the figure). To be more specific, the plunging molecular levels 
(dotted lines) and the BEC levels in the lowest gaseous adiabatic potential curve (solid lines) seem to repell each 
other at large a sc due to channel coupling. In summary, for our time-independent studies here, the strict adiabatic 
approximation describes the main characteristics of the energy level pattern as shown in Figs. || and |] with sufficient 
accuracy. To go beyond this static picture, however, the coupling between channels must be considered. 

Studies including channel coupling Jl7| , ^9| show that the primary recombination mechanism is related to the behavior 
of the adiabatic potential curves and the corresponding coupling elements at hyperradii R « 2 — 3a sc . Consider panel A 
of Fig. |]with a-sc = 2138a. u. and V} lo = 78kHz. For this situation, the criterion R w 2— 3a sc translates into 5.8— 8.8aho- 
At these hyperradii, the lowest lying gaseous Vt rap + U v curve reaches » 10 — 20?kj/j O , compared to the lowest gaseous 
energy level of w 5.3hu>ho- Thus, we speculate that we have entered a regime where recombination processes are 
largely suppressed due to the strong confinement. It would be interesting to probe this regime experimentally. 
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IV. COMPARISON WITH MEAN-FIELD TREATMENTS 



This section connects our microscopic studies with results obtained by mean-held treatments. Consider the HF 
treatment |11|, which results in an equation identical to the GP equation, except for a change from N to N — 1 in 
the interaction parameter. The HF treatment does not take advantage of the decoupling of the center-of-mass and 
the internal motion, and determines an approximate value for the total energy Ehf as well as the orbital energy or 
chemical potential. To compare our internal three-body energies Ej% (solid lines in Figs. | and|) with the HF energy 
Ehf, we subtract the exact energy associated with the lowest center-of-mass motion, E$ M = 1.5hu>h , from Ehf, 
and refer to this quantity as internal HF energy E l ^ F . The approximate HF treatment could either lead to an error 
in Eh F or in E^p 1 , or in both these quantities at the same time, and the subtraction of Eq from Ehf is therefore 
somewhat artificial |32|. Figures || and [|show the resulting internal HF energy, Effp, as dotted lines. 

For Vho = 780Hz, the internal HF energy E^ffp (dotted lines in Fig. ||) agrees favorably with the lowest gaseous-like 
internal three-body energy (solid line). As discussed in Sec. [II, Fig. [| includes only a small range of interaction 
parameters (N — \)a sc /a} lo . Thus, the good agreement between E l j£ and Ejgp is not surprising. 

For comparison, Fig. |^ shows the energy levels for a 100 times larger trapping frequency, namely Vho — 78kHz. 
This figure indicates two regions in which E l ^ F (dotted line) deviates from the lowest many-body energy level 
having gaseous character, the region around a cr and the large a sc region. For a sc < a cr , the HF equation does not 
have a solution. This mathematical fact is typically interpreted as "decay of the condensate" Jl]|l8|]. The three-body 
treatment reveals that the decay of the condensate is linked to the disappearance of a potential barrier (see, e.g., 
panels D and E of Fig. ||), which separates the molecular region at small R from the gaseous-like region at large R. 
When the potential barrier has disappeared, the states associated with this potential curve have purely molecular 
character. Thus, the energy levels in this potential curve change dramatically in the vicinity of a cr , whereas most of 
the energy levels belonging to other families change slowly around a sc . 

Our interpretation of the behavior of a BEC in the vicinity of a cr agrees with that suggested by Bohn et al. for an 
arbitrary number of particles pjjfl . Within the K- harmonic approximation, Bohn et al. [ |33| describe the interaction 
between each pair of particles through a shape-independent delta-function potential, and make some approximations 
about the behavior at small hyperradii, where the J-function potential leads to divergences. These approximations 
are confirmed to be well justified, based on our calculations with shape-dependent two-body interaction potentials, 
so that no additional approximations at small R need to be made. Taken together, these two studies highlight the 
usefulness of the hyperspherical radius coordinate in interpreting dynamics of a BEC (see also |34|]). 

Figure || also indicates discrepencies between E % g F and E l j£ at large |a sc |. For a sc — > oo, E l g F diverges and 
approaches unrealistically large values, whereas the three-body energies E™^ behave smoothly. Recall that the sign 
change of the scattering length from — oo to +00 (pole region) corresponds to the appearance of an additional two-body 
s-wave bound state. Coincidently, this corresponds to the existence of a hyperradial potential curve, which approaches 
the corresponding (negative) two-body binding energy in the absence of the trapping potential asymptotically. Around 
the pole region, this hyperradial potential curve is energetically close to the lowest potential curve with gaseous-like 
character, and therefore, coupling between these two potential curves with molecular-like and gaseous-like character, 
respectively, can occur. 

For positive a sc , the "standard" HF/GP treatment has been improved by including an additional mean-held 
term p |35| , |36[ | . This modified treatment leads to an improved description of BECs with small and medium inter- 
action parameters, however, exhibits divergences similar to that observed for the standard mean-field treatment as 
a sc — > 00. Therefore the modified mean-field treatment does not overcome the limitations of the standard mean-field 
treatment addressed above. 

To improve upon the GP (or equivalently, the HF) treatment, several groups fTzj-fhif developed mean- field formalisms 
within the shape-independent contact potential approximation. These approaches go beyond typically applied ap- 
proximations by introducing a molecular field, and are aimed at describing condensate physics in the presence of a 
Feshbach resonance. Application of these formalisms to calculating the energetics as a function of a sc , i.e., in the 
regions where a sc w a cr and |o ac | ~ > 00 would be useful. The present study may then be used to benchmark such 
treatments, and would thus provide a strong link between many-body and mean-field treatments. 

To make further contact with mean-field treatments, comparisons between our three-body excitation frequencies 
and those calculated within the random phase approximation can be made. This comparison will be discussed in a 
future publication. 
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V. CONCLUSIONS 



This paper discusses microscopic studies for three 87 Rb atoms interacting via a sum of simple two-body model 
potentials under external confinement. Our study is important in light of recent experiments of atomic gases trapped 
in an optical lattice with an occupancy of 1-3 atoms per lattice site . These experiments open the opportunity to 
systematically study strongly-interacting few-body systems, and therefore enter regimes discussed in the present paper 
from a theo retical microscopic point of view. The overall behavior of the zero temperature energy levels discussed 
in Sees. Ill and [V is independent of the exact shape of the two-body potential, though the detailed behavior of the 
energy levels and of the hyperradial potential curves at small hyperradii R depends on the shape of the two-body 
potential. The goal of this paper is to point out the gross behavior of the energy levels as a function of the two-body 
s-wave scattering length a sc . Therefore, we consider only one class of two-body po tent ials. Our main results are 
derived for a relatively large trapping frequency, Vho — 78kHz. As pointed out in Sec. Ill, dependencies on Vho exist, 
however, the gross features of the energy level scheme are independent of the magnitude of Vho- 

Our time-independent calculations, performed mostly within the adiabatic approximation, include molecular-like 
and gaseous-like states, whereas standard mean-field treatments only treat the latter states. Our study reveals a 



331]. Furthermore, the 
evels around the pole 



microscopic understanding of the decay of a BEC in agreement with a study by Bohn et al. 
present microscopic three-body study reveals a smooth non- diverging behavior of the energy 
region where |a sc | — > oo, and provides a detailed picture of the corresponding physics through analysis of hyperspherical 
potential curves. 

Section |v] compares our microscopic energy levels with those calculated at the HF level. While mean-field theories 
are commonly derived in the large particle limit, the HF equation is valid for any number of particles and density 
(though it might result in a more accurate description for systems with large N and low density). The energy level 
pattern for systems with more than three particles may differ in detail from our three-body energy levels presented 
here, however, the overall behavior will be similar and we expect our main conclusions to generalize to systems with 
more than N — 3 particles. 

Our study raises several discussion points. Most of all, we have only briefly discussed the effect of channel coupling 
on the behavior of the energy level scheme (see Fig. |^). To obtain detailed information about the dynamics of the 
condensate, inclusion of channel coupling, or at least a Landau-Zener-type analysis, is required. Such an analysis 
would allow the following question to be addressed (see Fig. ^): Is it possible to form a BEC with positive a BC that 
resides in a hyperradial curve other than the lowest curve with BEC character, and then decrease a sc until a sc < a cr 
(e.g. through the use of a Feshbach resonance) without destroying the condensate? Answers to this question certainly 
require a more detailed analysis, and will be left for future study. 

To extend the present work to many (more than four) particles is computationally infeasible. We hope, however, 
that our results can eventually be compared with modified mean-field treatments that implicitly include molecular 
states I^Q. Our study may then provide an additional link between many-body and mean-field physics. Lastly, our 
work relates to various many-body treatments using hyperspherical coordinates p3|j3^j45| , ^ | that neglect molecular 
physics. (An exception is the recent study by S0rensen et al. JIq|, which includes molecular states explicitly.) 
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FIG. 1. Upper panel: Adiabatic potential curves, U u , plus trapping potential, Vtrap ivho = 78kHz, aho = 731a.u.), as a 
function of the hyperradius R for Nt = 2 and a sc = 228a.u.. The dashed curves (y = 0, 1) describe molecular-like states while 
the solid curves describe gaseous-like states. Lower panel: Blowup of the gaseous region. In addition, dotted lines show the 
potential curves for vanishing interaction potential V. Note the logarithmic R scale. 



FIG. 2. Internal energy levels (solid lines) as a function of [arctan(a ac /ah )]/7r for three particles in a spherical trap 
with trapping frequency Vho = 78kHz (aho = 731a.u.). Dotted lines indicate the "internal" HF energy E l ^ F (see text). The 
"first cycle" (Nt = 1, see text) corresponds to [arctan(a sc /aho)]/7r £ [—0.5, —1.5], and the "second cycle" (iV;, = 2, see text) 
to [arctan(a sc /ah )]/7r 6 [—1.5,-2.5] ([arctan(a sc /ak )]/7r = — 0.5, — 1, —1.5, —2, —2.5 corresponds to |a sc | = oo, 0, oo, 0, oo). 
Labels A-F indicate systems with different scattering length for which Fig. ^| shows the corresponding adiabatic potential 
curves (A: a sc = 2138a.u.; B: a sc = 103a.u.; C: a sc = Oa.u.; D: a sc — — 31a.u.; E: a sc — — 130a.u.; and F: a sc = — 1970a.u.). 
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FIG. 3. Adiabatic potential curves, U„, plus trapping potential, V tra p, for three particles in a trap with v^o = 78kHz 
(ciho = 731a. u.) as a function of the hyperradius R. Solid lines indicate the numerical potential curves for A: a ac — 2138a.u.; 
B: cisc = 103a. u. (d-wave shape resonance); C: a sc = Oa.u.; D: a sc = — 31a.u. (first signature of decay); E: a sc = — 130a.u. 
(vanishing potential barrier); and F: a sc — — 1970a.u.. These scattering lengths are also indicated in Fig. H. The thick solid line 
indicates the lowest gaseous level (panel A-D), and the "decayed" condensate state (panel E and F), respectively. In addition, 
dotted lines show the potential curves for vanishing two-body interaction potential. 



FIG. 4. Internal energy levels (solid lines) as a function of [arctan(a ac /ah )]/7r for three particles in a spherical trap 
with trapping frequency Vho = 780Hz (flho = 7307a.u.) for Nt = 2. Diamonds indicate the internal energies at large negative 
a sc . A dotted line indicates the "internal" HF energy E^p (see text). 

FIG. 5. Internal adiabatic energy levels (dotted, solid and dashed lines) together with energies calculated using the 

coupled-adiabatic-channel approach including four channels (pluses) as a function of [arctan(a sc /a^ )]/7r for three particles 
in a spherical trap with trapping frequency «/; lo = 78kHz (flho = 731a.u.). Dotted lines indicate "plunging" molecular levels, 
while solid and dashed lines indicate energy levels in the lowest and in the second lowest gaseous adiabatic potential curve, 
respectively. 
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